{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# PyMOL Visualization\n",
    "PyMOL作为最流行的分子结构可视化的工具之一，在PyRosetta中内置了相关的服务可以直接和PyMOL进行工作联动。这使得我们可以实时观察结构和打分的变化。这大大有助于理解采样和计算过程中发了什么。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1 PyMOLRosettaServer的配置\n",
    "PyRosetta的一大特点在于可以使用PyMOL作为可视化的窗口，实时观测采样的过程。 首先安装PYMOL, 通过conda可以轻松获取。\n",
    "\n",
    "运行命令:\n",
    "```\n",
    "conda install -c schrodinger pymol\n",
    "```\n",
    "\n",
    "启动PyMOl命令:\n",
    "```\n",
    "pymol\n",
    "```\n",
    "\n",
    "图形界面菜单栏-点击-File>Edit pymolrc. 在其中添加:（路径根据自己PyRosetta安装的地址进行填写。）\n",
    "```\n",
    "run /opt/miniconda3/envs/pyrosetta/lib/python3.6/site-packages/pyrosetta/PyMOLRosettaServer.py\n",
    "```\n",
    "\n",
    "**重启PyMOL后如有以下字样说明配置成功**:\n",
    "\n",
    "PyMOL>run /opt/miniconda3/envs/pyrosetta/lib/python3.7/site-packages/pyrosetta/PyMOLRosettaServer.py\n",
    "\n",
    "PyMOL <---> PyRosetta link started!\n",
    "\n",
    "at 127.0.0.1 port 65000\n",
    "\n",
    "![title](./img/pymolrosettaserver.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1. 导入PDB生成Pose对象\n",
    "在前面章节中，我们已经学会如何得到Pose的4种方法，那就先来生成Pose，并发送到PyMOL中进行观察。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "PyRosetta-4 2020 [Rosetta PyRosetta4.conda.mac.python37.Release 2020.02+release.22ef835b4a2647af94fcd6421a85720f07eddf12 2020-01-05T17:31:56] retrieved from: http://www.pyrosetta.org\n",
      "(C) Copyright Rosetta Commons Member Institutions. Created in JHU by Sergey Lyskov and PyRosetta Team.\n",
      "\u001b[0mcore.init: {0} \u001b[0mChecking for fconfig files in pwd and ./rosetta/flags\n",
      "\u001b[0mcore.init: {0} \u001b[0mRosetta version: PyRosetta4.conda.mac.python37.Release r242 2020.02+release.22ef835b4a2 22ef835b4a2647af94fcd6421a85720f07eddf12 http://www.pyrosetta.org 2020-01-05T17:31:56\n",
      "\u001b[0mcore.init: {0} \u001b[0mcommand: PyRosetta -ex1 -ex2aro -database /opt/miniconda3/lib/python3.7/site-packages/pyrosetta/database\n",
      "\u001b[0mbasic.random.init_random_generator: {0} \u001b[0m'RNG device' seed mode, using '/dev/urandom', seed=-1937431081 seed_offset=0 real_seed=-1937431081 thread_index=0\n",
      "\u001b[0mbasic.random.init_random_generator: {0} \u001b[0mRandomGenerator:init: Normal mode, seed=-1937431081 RG_type=mt19937\n",
      "\u001b[0mcore.chemical.GlobalResidueTypeSet: {0} \u001b[0mFinished initializing fa_standard residue type set.  Created 980 residue types\n",
      "\u001b[0mcore.chemical.GlobalResidueTypeSet: {0} \u001b[0mTotal time to initialize 0.722618 seconds.\n",
      "\u001b[0mcore.import_pose.import_pose: {0} \u001b[0mFile './data/pose_demo.pdb' automatically determined to be of type PDB\n",
      "\u001b[0mcore.conformation.Conformation: {0} \u001b[0mFound disulfide between residues 2 11\n",
      "\u001b[0mcore.conformation.Conformation: {0} \u001b[0mcurrent variant for 2 CYS\n",
      "\u001b[0mcore.conformation.Conformation: {0} \u001b[0mcurrent variant for 11 CYS\n",
      "\u001b[0mcore.conformation.Conformation: {0} \u001b[0mcurrent variant for 2 CYD\n",
      "\u001b[0mcore.conformation.Conformation: {0} \u001b[0mcurrent variant for 11 CYD\n",
      "\u001b[0mcore.conformation.Conformation: {0} \u001b[0mFound disulfide between residues 5 25\n",
      "\u001b[0mcore.conformation.Conformation: {0} \u001b[0mcurrent variant for 5 CYS\n",
      "\u001b[0mcore.conformation.Conformation: {0} \u001b[0mcurrent variant for 25 CYS\n",
      "\u001b[0mcore.conformation.Conformation: {0} \u001b[0mcurrent variant for 5 CYD\n",
      "\u001b[0mcore.conformation.Conformation: {0} \u001b[0mcurrent variant for 25 CYD\n"
     ]
    }
   ],
   "source": [
    "from pyrosetta import pose_from_pdb,init,create_score_function\n",
    "init()\n",
    "pdb_pose = pose_from_pdb('./data/pose_demo.pdb')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 发送pose到pymol中\n",
    "from pyrosetta.rosetta.protocols.moves import PyMOLMover\n",
    "pymover = PyMOLMover()\n",
    "pymover.update_energy(True)\n",
    "pymover.apply(pdb_pose)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "此处截图可见，Pose已经发送至PyMOL窗口。\n",
    "![title](./img/sendpose.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2. 调用PyMOLMover进行能量分析\n",
    "PyMOLMover能量相关的api函数:\n",
    "- send_hbonds \n",
    "- send_energy/send_RAW_Energies\n",
    "- send_polars(注意:不兼容新版本的pymol)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\u001b[0mcore.scoring.etable: {0} \u001b[0mStarting energy table calculation\n",
      "\u001b[0mcore.scoring.etable: {0} \u001b[0msmooth_etable: changing atr/rep split to bottom of energy well\n",
      "\u001b[0mcore.scoring.etable: {0} \u001b[0msmooth_etable: spline smoothing lj etables (maxdis = 6)\n",
      "\u001b[0mcore.scoring.etable: {0} \u001b[0msmooth_etable: spline smoothing solvation etables (max_dis = 6)\n",
      "\u001b[0mcore.scoring.etable: {0} \u001b[0mFinished calculating energy tables.\n",
      "\u001b[0mbasic.io.database: {0} \u001b[0mDatabase file opened: scoring/score_functions/hbonds/ref2015_params/HBPoly1D.csv\n",
      "\u001b[0mbasic.io.database: {0} \u001b[0mDatabase file opened: scoring/score_functions/hbonds/ref2015_params/HBFadeIntervals.csv\n",
      "\u001b[0mbasic.io.database: {0} \u001b[0mDatabase file opened: scoring/score_functions/hbonds/ref2015_params/HBEval.csv\n",
      "\u001b[0mbasic.io.database: {0} \u001b[0mDatabase file opened: scoring/score_functions/hbonds/ref2015_params/DonStrength.csv\n",
      "\u001b[0mbasic.io.database: {0} \u001b[0mDatabase file opened: scoring/score_functions/hbonds/ref2015_params/AccStrength.csv\n",
      "\u001b[0mbasic.io.database: {0} \u001b[0mDatabase file opened: scoring/score_functions/rama/fd/all.ramaProb\n",
      "\u001b[0mbasic.io.database: {0} \u001b[0mDatabase file opened: scoring/score_functions/rama/fd/prepro.ramaProb\n",
      "\u001b[0mbasic.io.database: {0} \u001b[0mDatabase file opened: scoring/score_functions/omega/omega_ppdep.all.txt\n",
      "\u001b[0mbasic.io.database: {0} \u001b[0mDatabase file opened: scoring/score_functions/omega/omega_ppdep.gly.txt\n",
      "\u001b[0mbasic.io.database: {0} \u001b[0mDatabase file opened: scoring/score_functions/omega/omega_ppdep.pro.txt\n",
      "\u001b[0mbasic.io.database: {0} \u001b[0mDatabase file opened: scoring/score_functions/omega/omega_ppdep.valile.txt\n",
      "\u001b[0mbasic.io.database: {0} \u001b[0mDatabase file opened: scoring/score_functions/P_AA_pp/P_AA\n",
      "\u001b[0mbasic.io.database: {0} \u001b[0mDatabase file opened: scoring/score_functions/P_AA_pp/P_AA_n\n",
      "\u001b[0mcore.scoring.P_AA: {0} \u001b[0mshapovalov_lib::shap_p_aa_pp_smooth_level of 1( aka low_smooth ) got activated.\n",
      "\u001b[0mbasic.io.database: {0} \u001b[0mDatabase file opened: scoring/score_functions/P_AA_pp/shapovalov/10deg/kappa131/a20.prop\n",
      "\u001b[0mbasic.io.database: {0} \u001b[0mDatabase file opened: scoring/score_functions/elec_cp_reps.dat\n",
      "\u001b[0mcore.scoring.elec.util: {0} \u001b[0mRead 40 countpair representative atoms\n",
      "\u001b[0mcore.pack.dunbrack.RotamerLibrary: {0} \u001b[0mshapovalov_lib_fixes_enable option is true.\n",
      "\u001b[0mcore.pack.dunbrack.RotamerLibrary: {0} \u001b[0mshapovalov_lib::shap_dun10_smooth_level of 1( aka lowest_smooth ) got activated.\n",
      "\u001b[0mcore.pack.dunbrack.RotamerLibrary: {0} \u001b[0mBinary rotamer library selected: /opt/miniconda3/lib/python3.7/site-packages/pyrosetta/database/rotamer/shapovalov/StpDwn_0-0-0/Dunbrack10.lib.bin\n",
      "\u001b[0mcore.pack.dunbrack.RotamerLibrary: {0} \u001b[0mUsing Dunbrack library binary file '/opt/miniconda3/lib/python3.7/site-packages/pyrosetta/database/rotamer/shapovalov/StpDwn_0-0-0/Dunbrack10.lib.bin'.\n",
      "\u001b[0mcore.pack.dunbrack.RotamerLibrary: {0} \u001b[0mDunbrack 2010 library took 0.19348 seconds to load from binary\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "-43.15058616145932"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# 对pose先进行能量计算（必须执行）\n",
    "scoring = create_score_function('ref2015')\n",
    "scoring(pdb_pose)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 显示所有氢键的情况:\n",
    "pymover.send_hbonds(pdb_pose)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "![title](./img/pymol_hbnet.png)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 发送pose的能量着色到pymol中\n",
    "pymover.send_energy(pdb_pose)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "蓝色代表这个残基的能量越低(favorable energy)，残基着色越红代表能量越高(unfavorable energy)\n",
    "![title](./img/pymol_energy.png)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "除了整体着色以外，也可以支持仅对某个能量项进行着色，比如，我目前最关心分子内部的溶剂化能（看看哪些区域有很大暴露的疏水氨基酸，这会导致蛋白质折叠不稳定）。运行结果表示，β折叠片上的氨基酸侧链上的疏水原子有比较好的packing，暴露的疏水面积较α-螺旋上的少。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 发送某个能量打分项到pymol中\n",
    "pymover.send_energy(pdb_pose,'fa_sol')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "![title](./img/pymover_fasol.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 3. 调用PyMOLMover进行实时轨迹观察"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "除了静态的分析能量，在PyMOLMover中可以调用PyMOL Observer实时观察能量最小化过程\n",
    "\n",
    "动态历史记录相关的api函数:\n",
    "- keep_history\n",
    "- update_interval/update_energy\n",
    "\n",
    "**AddPyMOLObserver_to_conformation顾明思议，仅当Pose发生构象变化时，通过观察器提醒pymover更新并输出。\n",
    "\n",
    "**AddPyMOLObserver_to_energies，仅当能量发生变化时，通过观察器提醒pymover更新并输出。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "PyRosetta-4 2020 [Rosetta PyRosetta4.conda.mac.python37.Release 2020.02+release.22ef835b4a2647af94fcd6421a85720f07eddf12 2020-01-05T17:31:56] retrieved from: http://www.pyrosetta.org\n",
      "(C) Copyright Rosetta Commons Member Institutions. Created in JHU by Sergey Lyskov and PyRosetta Team.\n",
      "\u001b[0mcore.init: {0} \u001b[0mChecking for fconfig files in pwd and ./rosetta/flags\n",
      "\u001b[0mcore.init: {0} \u001b[0mRosetta version: PyRosetta4.conda.mac.python37.Release r242 2020.02+release.22ef835b4a2 22ef835b4a2647af94fcd6421a85720f07eddf12 http://www.pyrosetta.org 2020-01-05T17:31:56\n",
      "\u001b[0mcore.init: {0} \u001b[0mcommand: PyRosetta -ex1 -ex2aro -database /opt/miniconda3/lib/python3.7/site-packages/pyrosetta/database\n",
      "\u001b[0mbasic.random.init_random_generator: {0} \u001b[0m'RNG device' seed mode, using '/dev/urandom', seed=568641923 seed_offset=0 real_seed=568641923 thread_index=0\n",
      "\u001b[0mbasic.random.init_random_generator: {0} \u001b[0mRandomGenerator:init: Normal mode, seed=568641923 RG_type=mt19937\n",
      "\u001b[0mcore.import_pose.import_pose: {0} \u001b[0mFile './data/4jfx_peptide.pdb' automatically determined to be of type PDB\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "12.957104698417915"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# loading a pose;\n",
    "from pyrosetta import pose_from_pdb,init,create_score_function\n",
    "init()\n",
    "peptide_pose = pose_from_pdb('./data/4jfx_peptide.pdb')\n",
    "\n",
    "# 对pose先进行能量计算（必须执行）\n",
    "scoring = create_score_function('ref2015')\n",
    "scoring(peptide_pose)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "#AddPyMOLObserver(pose, True)\n",
    "from pyrosetta.rosetta.protocols.moves import AddPyMOLObserver_to_conformation\n",
    "from pyrosetta.rosetta.protocols.minimization_packing import MinMover\n",
    "from pyrosetta.rosetta.core.kinematics import MoveMap\n",
    "from pyrosetta.rosetta.protocols.moves import PyMOLMover\n",
    "pymover = PyMOLMover()\n",
    "pymover.send_energy(peptide_pose)\n",
    "pymover.apply(peptide_pose)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\u001b[0mcore.scoring.ScoreFunctionFactory: {0} \u001b[0mSCOREFUNCTION: \u001b[32mref2015\u001b[0m\n"
     ]
    }
   ],
   "source": [
    "# 创建pymol观察器;\n",
    "AddPyMOLObserver_to_conformation(peptide_pose, True)\n",
    "\n",
    "# 对多肽进行能量最小化;\n",
    "mp = MoveMap()\n",
    "mp.set_bb(True)\n",
    "mp.set_chi(True)\n",
    "min_mover = MinMover()\n",
    "min_mover.movemap(mp)\n",
    "min_mover.apply(peptide_pose)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "![Alt Text](./img/pymover_min.gif)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "# update energy;\n",
    "pymover.send_energy(peptide_pose)\n",
    "pymover.apply(peptide_pose)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "![Alt Text](./img/min.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "能量最小化之后，3号位的VAL氨基酸能量更低了（紫色）"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
